function [v] = LK_alg(I1,I2,lambda,mask,v_initial,num_iterations)
    v = v_initial;
    
    [Ix Iy] = ImageDerivatives(I1,I2);

    IxIy_sum = matrixSum(Ix.*Iy);
    
    A = [matrixSum(Ix.*Ix)+lambda,IxIy_sum ; IxIy_sum, matrixSum(Iy.*Iy)+lambda];
    
    for i=1:num_iterations
        %warp I2 towards I1 according current v
        [warpedI2 warpingMask] = warp(I2,v);
        
        %update mask
        warpingMask = warpingMask.*mask;
        
        %calculate new translation for estimating new velocity
        It = (warpedI2 - I1).*warpingMask;
        
        %add current solution to the accumulative velocity
        B = -[matrixSum(Ix.*It),matrixSum(Iy.*It)]';
        
        v = v + A\B;

    end
    
end
